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Abstract We present here a timing analysis of the fastest accreting millisec- 
ond pulsar IGR J00291+5934 using RXTE data taken during the outburst 
of December 2004. We corrected the arrival times of all the events for the 
orbital (Doppler) effects and performed a timing analysis of the resulting 
phase delays. In this way we find a clear parabolic trend of the pulse phase 
delays showing that the pulsar is spinning up as a consequence of accretion 
torques during the X-ray outburst. The accretion torque gives us for the 
first time an independent estimate of the mass accretion rate onto the neu- 
tron star, which can be compared with the observed X-ray luminosity. We 
also report a revised value of the spin period of the pulsar. 
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1 INTRODUCTION 



The so-called recycling scenario links two different classes of astronomical objects, namely 
the millisecond radio pulsars (usually found in binary systems) and the Low Mass X-ray 
* E-mail: burderi@mporzio.astro.it 



2 L. Burderi et al. 

Binaries (hereafter LMXBs), or, at least, a subgroup of them. The leading idea of this 
scenario is the recycling process itself, during which an old, weakly magnetized, slowly 
spinning neutron star is accelerated by the accretion of matter and angular momentum 
from a (Keplerian) accretion disk down to spin periods in the millisecond range. In 
this way, at the end of the accretion phase, the neutron star rotates so fast that it is 
resurrected from the radio pulsar graveyard, allowing the radio pulsar phenomenon to 
occur again despite the weakness of the magnetic field. 

Although this scenario was firstly proposed long time ago (see e.g. Bhattacharya & 
van den Heuvel 1991 for a review), the most cmbarassing problem was the absence of 
coherent pulsations in LMXBs. Only recently, the long seeked for millisecond coherent os- 
cillations in LMXBs have been found, thanks to the capabilities (the right combination of 
high temporal resolution and large collecting area) of the RXTE satellite. In April 1998, 
a transient LMXB, SAX J1808. 4-3658, was discovered to harbour a millisecond pulsar 
(P sp i n ~ 2.5 ms) in a compact (P rb — 2 h) binary system (Wijnands & van der Klis 
1998; Chakrabarty & Morgan 1998). We now know seven accreting millisecond pulsars 
(Wijnands 2005; Morgan et al. 2005); all of them are X-ray transients in very compact 
systems (orbital period between 40 min and 4 h), the fastest of which (P sp m — 1.7 ms), 
IGR J00291+5934, has been discovered in December 2004 (Galloway et al. 2005, here- 
after G05). Timing techniques applied to data of various accreting millisecond pulsars, 
spanning the first few days of their outbursts, allowed an accurate determination of their 
main orbital parameters. However, only a few attemps have been made to determine the 
spin period derivative (Chakrabarty et al. 2003; Galloway et al. 2002). 

In this paper we apply an accurate timing technique to the fastest currently known 
accreting millisecond pulsar, IGR J00291+5934, in the hope of constraining the pre- 
dictions of different torque models with good quality experimental data. Our results 
indicate quite clearly that a net spin up occurs during the December 2004 outburst of 
IGR J00291+5934 (see also Falanga et al. 2005) and that the derived torque is in good 
agreement with that expected from matter accreting from a Keplerian disk. 

2 THE TIMING TECHNIQUE 

In standard timing techniques (see e.g. Blandford & Teukolsky 1976) the predicted arrival 
time of a given pulse is computed using a first guess of the parameters of the system, and 
the difference between the experimental and predicted arrival times, namely the residuals, 
are fitted with a linear multiple regression of the differential corrections to the parameters. 
This means that the differential correction to orbital parameters, source position in the 
sky, spin frequency and its derivative, are computed simultaneously. This technique has 
the obvious advantage to give a self-consistent solution, where all the correlations in the 
covariance matrix of the system are fully taken into account. However, the convergency of 
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the fit is not always guaranteed and - especially on long temporal baselines - convergence 
to secondary minima could occur. 

On the other hand, if the orbital period is much shorter than the timescale on which 
the spin period derivative is expected to produce a significant effect, we can demonstrate 
that a different timing technique is more effective in determining the spin period deriva- 
tive. This technique relies on the fact that the delays in the arrival times produced by 
the orbital corrections are effectively decoupled from those caused by the spin evolution. 
The technique proceedes as follows: in order to obtain the emission times, t eia , the arrival 
times of all the events, i arr , are firstly reported to the Solar system barycenter adopting 
the best estimate of the source position in the sky, then corrected for the delays caused by 
the binary motion using the best estimate of the orbital parameters through the formula: 

" 2tt 



^cm — ^arr SID. 



(tan - T*) 



(1) 



_ ^orb 

where x = asinz/c is the projected semimajor axis in light seconds, and T* is the 
time of ascending node passage at the begining of the observation. In the following, for 
simplicity, we use t instead of t em . The differential of this expression with respect to the 
orbital parameters allows to calculate the uncertainties in the phases, cr^ or b, caused by 
the uncertainties in the estimates of the orbital parameters. In a similar way, we have 
computed the uncertainties in the phase delays, a^pos, caused by the uncertainties on 
the estimates of the source position in the sky. In this case we can estimate v fitting the 
measured phase variations, while the uncertainties in the adopted orbital parameters and 
source position will result in a "timing noise" of amplitude fj0 par = (c| or b + ^pos) 1 ^ 2 
(see Burderi et al. 2005 for details). 



3 OBSERVATIONS AND DATA ANALYSIS 

IGR J00291+5934 was observed by RXTE between 2004 December 3 and 21. In this paper 
we report on the data between December 7 and 21 taken from a public ToO. We mainly 
use data from the PCA for timing anaysis and data from PCA and HEXTE for spectral 
analysis. The arrival times of all the events were converted to barycentric dynamical 
times at the solar system barycenter. The position adopted for the source was that of 
the proposed radio counterpart (which is compatible with that of the proposed optical 
counterpart, see Rupen et al. 2004). We corrected the arrival times of all the events for 
the delays caused by the binary motion using eq. with the orbital parameters given 
in G05. 

Adopting the uncertainties in the estimates of the orbital parameters given in G05, 
the positional uncertainty of 0.04" radius reported by Rupen et al. (2004) we obtain: 
c^orb 0.01, u pos Sa 0.006 where we have maximized sin and cos functions with 1, and 
used t — To Ss 7 days. Therefore, we expect that the uncertainties in the orbital parameters 
and source position will cause a "timing noise" not greater than cr0 par x P sp i n ~ 0.02 ms. 
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To compute phases of good statistical significance we epoch folded each interval of 
data in which the pulsation was significantly detected at the spin period given in G05 with 
respect to the same reference epoch, To, corresponding to the bcgining of our observations. 
The fractional part of the phase was obtained fitting each pulse profile with a sinusoid of 
fixed period. To compute the associated errors we combined the statistical errors derived 
from the fit, cr^ sta t, with the errors <70 par . 

In order to derive the differential correction to the spin frequency, A^o, and its deriva- 
tive, j>o, at the time To we have to derive a functional form for the time dependence of the 
phase delays. We started from the following simple assumptions: i) the bolometric lumi- 
nosity L is a good tracer of the mass accretion rate M via the relation L = ((GM/R)M, 
where £ < 1, and G, M, and R are the gravitational constant and the neutron star mass 
and radius, respectively, ii) The matter accretes through a Keplerian disk truncated at 
the magnetospheric radius, i? m oc M~ a , by its interaction with the (dipolar) magnetic 
field of the neutron star. At R m the matter is forced to corotate with the magnetic field 
of the neutron star and is funneled (at least in part) towards the rotating magnetic 
poles, thus causing the pulsed emission. For standard disk accretion a = 2/7. Indeed 
we considered two extreme cases, namely a = 2/7 and a = 0, since a location of i? m 
independent of M has been proposed (see, e.g., Rappaport, Fregeau, and Spruit, 2004). 
iii) The matter accretes onto the neutron star its specific Keplerian angular momentum 
at i? m , I — (GMi? m ) 1 ' 2 , thus causing a material torque — I x M. A firm upper limit 
to this torque is given by the condition < £ max x M, with £ max — {GMRqo) 1 / 2 , 
where Rco — 1-50 x 10 8 m 1/<3 j/~ 2 / 3 is the corotation radius (namely the radius at which 
the Keplerian frequency equals the spin frequency v of the neutron star and beyond 
which accretion is centrifugally inhibited), and m — M/M Q . We do not consider any 
form of threading of the accretion disk by the magnetic field of the neutron star (see e.g. 
Rappaport, Fregeau, and Spruit 2004 for a description of the magnetic threading), which 
implies that the only torque acting during accretion is r^-. 

In these hypothesis the spin frequency derivative is v = £M/(2ttI), where I is the 
moment of inertia of the neutron star and we have neglected any variation of / caused 
by accretion. If M = M(t), we have v{t) = (2tt /)- 1 4M (M(i)/M ) 1 " a/2 , where £ = 
(GMi? m0 ) 1//2 , and i? m0 and Mo are i? m and M at t = T , respectively. For the a = 
case we assumed £o — ^ max . 

Since we assumed M(t) oc L(t), to determine the temporal dependence of M(t) we 
studied the energy spectra of the source for each continuous interval of data combining 
PCA and HEXTE data. All the spectra are well fitted with a model consisting of a 
power law with an exponential cutoff plus thermal emission from a Keplerian accretion 
disk modified by photoelectric absorption and a Gaussian iron line. In order to derive 
L(t) for each spectrum we made the simple assumption L(t) oc F/ 3 _i 50 \(t), which is the 
unabsorbed flux in the RXTE PCA plus HEXTE energy band (3-150 keV). A good fit of 
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Table 1 Orbital and spin parameters of IGR J00291+5934. 
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Parameter 



G05 



This work 



Projected semimajor axis, a\ sin i (lt-ms) 
Orbital period, P or b (s) 

Epoch of ascending node passage, T* (MJD) 

Eccentricity, e 

Spin frequency, i>o (Hz) 

Spin frequency derivative, z>o (H z / S ) 

Spin frequency derivative, i>o (Hz/s) (a = 0) 

Spin frequency derivative, i>o (Hz/s) (a = 2/7) 

Spin frequency derivative, i>o (Hz/s) (v = constant) 

Epoch of the spin period, To (MJD) 



64.993(2) 

8844.092(6) 

53345.1619258(4) 

< 2 x 1(T 4 (3 a) 
598.89213064(1) 

< 8 x 10" 13 (3 a) 



598.89213053(2) 

1.17(0.16) x 10" 12 
1.11(0.16) x 10" 12 
0.85(0.11) x 10" 12 
53346.184635 



Errors are given at la confidence level. 



Ft 3 _ 150 }(t) vs t between December 7 and 14 (Ai b s ~ 7.3 days) is given by the expression 
■£(3-150) (*) = -£(3-150) X [1 - (* - T )/t B ] with t B = 8.4 ±0.1 days, where F( 3 _ 150 ) is the 
unabsorbed flux at t = To- Therefore we have v(t) — vq x [1 — (t — To)/t B ] 1 ~ a / 2 , where 
the spin frequency derivative at t = To is vo = (27r/) _1 (o Mq. We have therefore fitted 
the phase delays with the function: 



b - Ai/ it - T ) - -i> (i - T ) 2 x 



G 



(t - Tp) 
t B 



(2) 



Using the best fit value for Ai'o we computed the improved spin frequency estimate and 
repeated the same procedure described at the begining of this paragraph, folding at the 
new estimate of the spin period. The new phases were fitted with eq. J2J). In this case, 
Aj/q was fully compatible with zero. These phases are plotted versus time in Figure ^ 
(upper panel) together with the residuals in units of a with respect to eq. (lower 
panel). The best fit estimates of vq and Vo are reported in Table 1 for three values of 
a, namely a = which correspond to a location of R m independent of the accretion 
rate (cfr. the model of Rappaport Fregeau, and Spruit 2004 in which R m = Rco f° r 
any M), the standard case a = 2/7 which corresponds to R m proportional to the Alfven 
radius, and a = 2 which has been given for comparison purposes and corresponds to a 
parabolic trend, expected in the case of constant M. The statistics is not good enough 
to distinguish between these three possibilities. 

From the best-fit value of the spin frequency derivative vq we can compute the mass 
accretion rate at t = To through the formula: M_io = 5.9 x !>_ 13 I i5 m~ 2 / 3 (i?co /Smo) 1 ' 2 , 



where M_io is Mq in units of 10 10 M Q yr 1 , J>_i3 is j> i n units of 10 13 s 2 , and 



15 



6 



L. Burderi et al. 



0.2 
0.15 
0.1 

t 

0.05 


2 

1° 
•2 

■4 

Time(MJD-53346.18) 

Fig. 1 Pulse phases computed folding at the spin period reported in Table 1 
and plotted versus time together with the best fit curve for a — 2/7 shown as a 
solid line (upper panel) and residuals in units of a (lower panel). Note that the 
linear term is fully compatible with zero. 



is I in units of 10 45 g cm 2 . In the following we will adopt the FPS equation of state for 
the neutron star matter for m — 1.4 and the spin frequency of IGR J00291+5934 which 
gives J 45 = 1.29 and R = 1.14 x 10 6 cm (see e.g. Cook, Shapiro & Teukolsky 1994). In 
order to compare the experimental estimate of Mo with the observed X-ray luminosity, 
we have to derive the bolometric luminosity L(t) from the observed flux F( 3 _ 150 ){t). To 
this end we consider the spectral shape at t = To in more detail. 

The power law is the dominant spectral component. This component presumably 
originates in an atmosphere of small optical depth just above each polar cap (see e.g. 
Poutanen & Gierlinski 2003; Gierlinski & Poutancn 2005), thus we neglect, to first order, 
any effect of the inclination of the emitting region with respect to the observer. On 
the other hand, we observe a single-peaked pulse profile, which means that we only see 
the emission from one of the polar caps (e.g. Kulkarni & Romanova 2005). We have 
therefore multiplied by 2 the unabsorbed flux of the power law in order to take into 
account the emission of the unseen polar cap. Assuming isotropic emission, we computed 
Lpl — 2Fpl (o.oo) x 47tgP = 1.5^0 3 x 10 37 d^ kpc erg/s, where ciskpc is the source distance 
in units of 5 kpc. The uncertainty on the luminosity has been evaluated propagating 
the uncertainties on the spectral parameters threated as they were independent on each 
other. 

The second component is the thermal emission from a Shakura-Sunyaev accretion 
disk. The fraction of the total luminosity that is emitted by the disc is given by the ratio: 
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0.5R/R m0 = 0.34. In this hypothesis L BB = 0.39/(1 - 0.39) x L PL0 = 9.6 x 10 36 d 2 5kpc 
erg/s. The total bolometric luminosity is therefore Lq — 2.46io'29 x 10 37 ^5 kpc erg/s. If we 
compare this luminosity with the mass accretion rate inferred from the timing analysis, 
we obtain an estimate of the source distance, which is: d = 9.47^ 2 'i kpc. 

4 CONCLUSIONS 

We have analysed RXTE data of the fastest known accreting millisecond pulsar, IGR 
J00291+5934, during the period 7—14 December, 2004. We report a revised estimate of 
the spin period and the spin period derivative. The source shows a strong spin-up, which 
indicates a mass accretion rate of about 8.5 x 10~ 9 M Q yr _1 . Comparing the bolometric 
luminosity of the source as derived from the X-ray spectrum with the mass accretion rate 
of the source as derived from the timing, we find a good agreement if we place the source 
at a quite large distance between 7 and 10 kpc. Note that 10 kpc is close to the outer 
edge of our Galaxy in the direction of IGR J00291+5934. Another possibility is that part 
of the luminosity of the system is not observed because emitted in other energy bands 
or because of occultation effects (which may be favoured if indeed the source is highly 
inclined) . 
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